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ABSTRACT 

We present the results of high precision, high resolution (R^68000) optical observations of the 
short-period (4d) eccentric binary system Alpha Virginis (Spica) showing the photospheric line-profile 
variability that in this system can be attributed to non-radial pulsations driven by tidal effects. 
Although scant in orbital phase coverage, the data provide S/N>2000 line profiles at full spectral 
resolution in the wavelength range AA400CE8500A, allowing a detailed study of the night-to-night 
variability as well as changes that occur on ^2 hr timescale. Using an ab initio theoretical calculation, 
we show that the line-profile variability can arise as a natural consequence of surface flows that are 
induced by the tidal interaction. 

Subject headings: stars: individual (a Vir, Spica) — stars: binaries (including multiple): close - 
line: profiles — stars: atmospheres — stars: rotation — techniques: spectroscopic 



1. INTRODUCTION 

a Virginis (Spica, HD 116658) is a double- lined spec- 
troscopic binary consisting of two early B-type stars in a 
short-period (~4 d) eccentric orbit. The Bl III-IV pri- 
mary was discovered to be a (3 Cephei-type star in early 
observations, with a pulsation period of 0.1738 d in both 
the light curve (Shobbrook et al. 1969) and radial ve- 
locity (RV) variations (Smak 1970). Although the light 
variations became undetectable a short time after they 
were first detected (Lomb 1978; Sterken et al. 1986), the 
line-profile variability persists (Smith 1985a, 1985b; Rid- 
dle 2000). Its most striking characteristic is the presence 
of discrete absorption and emission features travelling 
from the blue towards the red wing of the absorption 
line. 

The first detection of this type of periodic line-profile 
variability was made by Walker et al. (1979) in the 
rotationally-broadened 09.5V star £ Ophiuchi. In their 
seminal study, Vogt & Penrod (1983) showed that this 
variability could be explained by non-radial pulsations 
(NRP). Subsequent observations show a prevalence of 
this type of line-profile variability (c.f. Baade 1984; Gies 
& Kullavanijaya 1988; Reid et al. 1993; Fullerton et al. 
1996; Rivinius et al. 2001; Uytterhoeven et al. 2001). 
The numerous models based on the NRP mechanism are 
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reviewed by Townsend (1997a, 1997b), some of which in- 
corporate the effects of rapid rotation. A general feature 
of these models is to assume a prescription for the sur- 
face velocity field of the pulsating star, from which the 
perturbed photospheric absorption line-profiles are com- 
puted. 

A very different approach was used by Moreno et al. 
(2005), which calculates the surface velocity field from 
first principles. This calculated (but not assumed) field 
is then projected along the line-of-sight to the observer to 
produce photospheric absorption line profiles. The one 
basic approximation in this model is to assume that only 
the external layer of the star oscillates in response to the 
forces in the system, while the interior region is assumed 
to rotate rigidly. This one-layer approximation is essen- 
tially equivalent to the assumption that the surface layer 
behavior is primarily controlled by the binary companion 
perturbations. 

The advantages of this method are threefold: 1) it 
makes no a priori assumption regarding the mathemati- 
cal formulation of the tidal flow structure since the veloc- 
ity field v is derived from first principles; 2) the method is 
not limited to slow stellar rotation rates nor to small or- 
bital eccentricities; and 3) it is computationally inexpen- 
sive. It is not clear to what extent the response of the in- 
ner layers may affect the line-profile calculation because 
the amplitudes of induced oscillations decrease sharply in 
deeper layers (Dolgivon & Smel'Chakova 1992). In this 
paper we confront the results obtained from the Moreno 
et al. (2005) method with very high-quality observa- 
tional data of Spica. We show that the azimuthal velocity 
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Fig. 1. — The top panel shows the signal-to-noise measurements output by the Libre-Esprit package as a function of wavelength after 
averaging all the individual polarimetric exposures. The bottom panel shows the night-to-night differences in the signal-to-noise for order 
44 (the 18th) at 5150A. 



components of the tidal perturbations (the "tidal flows" ) 
play a dominant role in determining the line-profile vari- 
ability as produced using this theoretical framework. 

The structure of this paper is as follows: In Sections 2 
and 3 we describe the observations and the characteris- 
tics of Spica's spectrum and variability; in Section 4 we 
present the results of the model calculations of Spica's 
surface velocity field and line-profile variability; the con- 
clusions are presented in Section 5; and details of the 
data reduction process and a discussion of the local line 
profile effects and size of the perturbations on the stellar 
surface are given in the Appendices. 



TABLE 1 
ESPaDOnS Observations 



Date & 


N 


S/N 


S/N 


S/N 


JD 


Orbital 


Time UT 




4270A 4720A 6660A 




Phase 


20 12:32:28 


8 


1358 


1685 


1580 


46.0225 


0.3776 


21 14:23:35 


8 


1592 


1978 


1868 


47.1000 


0.6523 


22 12:44:32 


8 


1660 


2039 


1892 


48.0309 


0.8838 


23 13:17:38 


8 


1401 


1753 


1662 


49.0539 


0.1387 


24 14:21:52 


8 


1670 


2074 


1933 


50.0985 


0.3987 


25 14:45:23 


8 


1360 


1723 


1633 


51.1148 


0.6516 


26 09:34:45 


8 


1493 


1856 


1747 


51.8991 


0.8472 


27 11:30:58 


16 


2349 


2897 


2659 


52.9798 


0.1167 


28 11:16:02 


16 


1933 


2386 


2202 


53.9695 


0.3623 



This Table shows the UT day and time of the March 2008 obser- 
vations, the number of exposures summed, the flux based estimate 
of the signal-to-noise ratio output by the Libre-Esprit package for 
4270, 4720 and 6660A, the Julian date - 2454500 and the orbital 
phase of the binary at the time of observation. 



2. OBSERVATIONS AND DATA REDUCTION 

The observations were performed in queue mode at the 
Canada France Hawaii 3.6m telescope (CFHT) with the 
ESPaDOnS spectropolarimeter. ESPaDOnS is a fiber- 
fed cross-dispersed echelle covering 3700A to 10480A 111 
a single exposure with 40 spectral orders on the EEV1 



detector with a nominal spectral resolution of R=68,000 
in polarimetric mode. The instrument has a dedicated 
reduction package, Libre-Esprit, that automatically pro- 
cesses the data. This script does typical flat fielding, 
bias subtraction, wavelength calibration using both cali- 
bration lamp and Fabry-Perot frames and optimal spec- 
tral extraction (see Donati et al. 1997). Our integration 
times were 6 seconds per exposure in polarimetric mode. 
This mode uses a Wollaston prism to produce two sep- 
arate orthogonally polarized spectral orders on the de- 
tector. A polarimetric sequence consists of two groups 
of 4 individual exposures at changing wave-plate an- 
gles. Each group of four exposures is processed by Libre- 
Esprit into average intensity and polarization spectra. 
The doubling of the spectral orders when combined with 
the intrinsically broadened spatial point-spread-function 
caused by the image sheer gives a peak signal-to-noise ra- 
tio above 1500 in individual spectral pixels without sat- 
uration. Thus, a polarimetric sequence consisting of 8 
exposures (16 individual spectra) spatially broadened by 
the image-sheer gives very high precision at full spectral 
resolution. 

On each night during 2008 March 20-March 26 we ob- 
tained one full polarimetric sequence. On March 27 and 
28 two complete sequences were obtained. Table [T] lists 
the date and time of the exposures, the number of im- 
ages combined to produce the spectrum, the signal to 
noise of for various lines in this combined spectrum, the 
Julian date and finally the corresponding orbital phase. 
The orbital phase was computed using P and To listed 
in Table [2] from Hcrbison-Evans et al. (1971). A phase 
of zero corresponds to periastron passage. 

ESPaDOnS is optimized for a broad wavelength cov- 
erage, having a high efficiency from roug hly 4000-8000A. 
Libre-Esprit automatically outputs an estimate of the 
signal-to-noise for each spectral order based on the av- 
erage flux in each order. Figure [T] shows the wavelength 
and temporal dependence of the signal-to-noise ratio for 
our data set. The peak signal-to-noise occurs in order 
17 around 5150A. The ratio is above half the peak value 
over the 4000-8500A range. Table [T] shows that the pre- 
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cision obtained for all the lines presented in this paper 
are similar. The ratio only changes by 35% from 4270 to 
6660A. 
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Fig. 2. — Amplified portions of the March 26 (top) and 28 (bot- 
tom) spectra showing the relative strengths and locations of the 
absorptions from mi and rrt2, as well as the principal identifica- 
tions of each line. The tick marks coincide with the laboratory 
wavelength of the identification. 



All cross-dispersed echelle spectrographs have multiple 
spectral orders with wavelength-dependent order over- 
lap. Typically shorter wavelength orders overlap signif- 
icantly while longer wavelength orders leave order gaps. 
The spectral orders must be individually continuum- 
normalized before being spectrally merged to account 
for the wavelength dependent sensitivity. The Libre- 
Esprit package produces both continuum-normalized and 
un-normalized spectra. In this package, each individual 
spectral order is rectified using a high-order polynomial 
response function to the continuum regions. However, 
this polynomial fit was not accurate enough for this work. 
The continuum normalization was done on each un- 
normalizcd individual exposure after orders were spec- 
trally merged. This merging and continuum mis-fitting 
can be problematic for accurate spectroscopic analysis. 
The analysis presented in this paper avoids lines falling 
near order boundaries and in regions where the order 
overlap is incomplete. This continuum fitting provides 
sufficient accuracy for this work (see appendix and ES- 
PaDOnS instrument description on the web). 

We have created IDL scripts to deal with the non- 
uniform overlapping spectral regions by binning the in- 
tensity values into essentially constant wavelength inter- 
vals of AA =0.lA. This wavelength interval is typically 
4 times the individual spectral-pixel size. These scripts 
result in a monotonic and nearly-regular wavelength cov- 
erage as well as increased signal-to-noise per wavelength 
bin resulting from the 2-12 pixel averaging. The proper- 
ties of the ESPaDOnS spectral coverage and a detailed 
description of the effects of this script are contained in 
the appendix. As an example, Figure [2] shows a plot of 



the March 20 merged and rectified spectrum where the 
very high signal-to-noise profiles can be seen. 



TABLE 2 

Spica Parameters 



Parameter 


Set "HE" 


Set "Aufden" 


mi (M Q ) 


io.gio.9 1 


10.25±0.68 4 


m 2 (Mq) 


e.sio.7 1 


6.97±0.46 4 ' a 


Ri (R©) 


s.iio.5 1 


7.40±0.57 4 ' a 


R 2 (R©) 




3.64±0.28 4 


P (days) 


4.014597 1 


4.0145898 


To (JD) 


2440678.09 1 


2440678.09 


e 


0.146 1 


0.067±0.014 3 


i (deg) 


66±2 1 


54±6 4 


w(deg) at T 


138 ±15 


140 ±10 4 


Aps. Period (yrs) 


mill 1 


135±15 4 


visini (km/s) 


161±2 2 


161±2 2 


v ro ti (km/s) 


176±5 


199±5 


visini (km/s) 


70±5 2 


70±5 2 


v rot 2 (km/s) 


77±6 


87±6 


#>(ml) 


1.3 


1.88±0.19 


#>(m2) 




1.67±0.5 



This Table lists the parameters of the Spica binary system from 
Herbison-Evans et al. 1971 (H-E 1971) and Aufdenberg et al. 2009 
(Aufden 2009) Reference 1 is Herbison-Evans et al. 1971, reference 2 
is Smith, 1985a, reference 3 is Riddle, R.L. 2000; reference 4 is Auf- 
denberg et al. 2009; note a these are the radii at the poles. 



3. SPECTRUM AND LINE-PROFILE VARIABILITY 

The spectrum of Spica is very similar to that described 
by Riddle (2000). It contains two sets of photospheric 
absorptions lines, each typical of an early B-type spec- 
trum. The components arising in the primary star, mi 
are significantly stronger than the lines arising in the sec- 
ondary component, m^. Figure [2] displays a fragment of 
the spectrum at two orbital phases close to elongations, 
where the two sets of lines are clearly resolved. The 
components are indicated in this Figure, along with the 
identifications of the strongest absorption lines. 

The most noteworthy features of Spica's spectra are 
the following: 1) The width of mi's lines is consistent 
with the value of vsini =161 km/s reported by Smith 
(1985), while that of m,2 is, indeed, significantly smaller. 
2) We are able to clearly identify He II 4685.7. Although 
its blue wing is blended with Si III 4683.80, the relatively 
small shift (—20 km/s) of the centroid of the feature with 
respect to mi 's reference frame indicates that He II is the 
dominant contributor to this blend. If both Si III and 
He II contributed equally, the expected shift would be 
—64 km/s. We also confirm the absence of He II 4541.6, 
thus concurring with Riddle (2000) and Aufdenberg et al. 
(2009) that mi is a B0. 5 star. 3) All weak lines present 
in mi's spectrum display the discrete narrow features, 
while lines having larger opacities, such as the H-Balmer 
series and the He I lines do not. We do not measure 
any significant difference in the equivalent widths of m.2 
lines at the two opposite elongation phases illustrated in 
Figure H 

3.1. Radial Velocities 
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Fig. 3. — Line centroids measured at different depths on the 
Si III 4552.62 line. The observation date and orbital phase are 
listed in the Figure. Each vertical dash corresponds to the ccntroid 
measured at that depth. The vertical dotted line is the average of 
all centroids shown in this Figure. The horizontal line shows the 
range of these calculated centroids. Any individual measurement 
can differ from the average by up to 16 km/s. This illustrates 
the source of uncertainty in the RV measurements. Note that the 
spectrum of the secondary is clearly separated yet the intrinsic line 
profile variability of the primary produces a very large systematic 
uncertainty. 



The first step in the line-profile analysis consists of 
transforming the wavelength scale to a velocity scale in 
the frame of reference of mi, which requires knowledge 
of its orbital motion. There are currently at least two 
determinations of the stellar and orbital parameters of 
Spica, the first from Herbison-Evans et al. (1971), the 
second from Aufdenberg et al (2009). The parameters we 
used are listed in Table[2]in two columns, the first derived 
from the results from Herbison-Evans et al. (labeled "H- 
E 1971" and the second derived from those of Aufdenberg 
et al. (labeled "Aufden 2009"). 

For the analysis presented in this paper, we chose the 
lines of Si III 4552.62 and He II 4685.7 (+Si III) for the 
following reasons: 1) they are relatively isolated from 
other neighboring lines; 2) they are relatively free from 
contamination by other atomic transitions; and 3) the 
contribution from 7712 is clearly visible in Si III, but ab- 
sent from He II. We also note that the equivalent widths 
of these Si III lines are insensitive to effective tempera- 
ture. 

The strong line-profile variability in Spica requires that 
the measurement of spectral lines be performed with 
great caution. We exemplify the problem in Figure [3] 
with the Si III 4552. 62A line observed on 26 March. 
The location of the line centroid is indicated for mea- 
surements performed at different intensity levels; as the 
measured intensity level approaches the continuum, the 
centroid of the line becomes increasingly more negative. 
The RV measured near the minimum of the line differs by 
32 km/s from the RV measured at the continuum level. 
This particular example is very illustrative because the 
contribution from m-2, clearly visible in the plot, cannot 



TABLE 3 
Radial Velocities 



Day 


JD 


Phase 


Si III 




He II 




He I 


m 2 


Na I 


Mar 






RV 


err 


RV 


err 


RV 


EW 


ISM 


20 


546.022 


0.377 


94 


3 


64 


4 


-115 


0.086 


-7 


21 


547.100 


0.652 


-77 


11 


-91 


1 






-8 


22 


548.031 


0.878 


-87 


11 


-118 


10 


179 


0.11 


-8 


23 


549.054 


0.132 


73 


6 


44 


3 


-135 


0.12 


-9 


24 


550.098 


0.393 


84 


3 


43 


2 






-10 


25 


551.115 


0.646 


-79 


9 


-96 


1 






-8 


26 


551.899 


0.841 


-106 


16 


-141 


8 


193 


0.10 


-10 


27 


552.980 


0.110 


68 


3 


39 


i 






-9 


28 


553.970 


0.357 


103 


2 


71 


2 


-139 


0.11 


-8 



This Table shows the radial velocities measured for each obser- 
vation with respect to laboratory wavelengths. The columns show 
the March 2008 observation date, Julian date (JD-2454000), orbital 
phase of observation, radial velocity (RV) and uncertainty (err) of 
the primary star for Si III (4552. 2A) and He II (4685. 7A) lines, the 
secondary companion He I (6678. 15A) radial velocity and equiv- 
alent width in angstroms and finally the Na I D ISM line radial 
velocity. 



contribute in any manner to the shape of mi's profile. 
Thus, it is very difficult at this stage to decide what por- 
tion of the line profile most reliably represents orbital 
motion and this situation introduces a high degree of un- 
certainty in the determination of the stellar and orbital 
parameters. 

Given this uncertainty, we measured Si III and He II 
4685. 7A at as many intensity levels as possible, without 
reaching levels where m,2's absorption is evident. The 
average value of these measurements is listed in Table El 
together with an uncertainty that corresponds to 1/2 of 
the total range in RV values measured for the line. The 
resulting RVs are plotted in Figure|U Also plotted in this 
Figure are the RV curves corresponding to the Aufden- 
berg et al. (2009) parameters (solid) for three different 
values of u>, the argument of periastron, which defines 
the location of line of apsides with respect to the plane 
of the sky. Aufdenberg et al. (2009) give to =140°±10° 
at T = JD 2440678.09, as well as the apsidal period 
Paps =135±15 years. Thus, for the current observations, 
W2008 lies between 248° and 257°. The measured radial 
velocities indicate that W2008 =255°, as illustrated in Fig- 
ure [H Since this lies within the range of acceptable val- 
ues, we adopt W2008 =255° for this study. 

3.2. Line-profile variability 

Line profile variability in Spica was first analyzed in 
detail by Smith (1985a, 1985b) in a closely-spaced set 
of data obtained over the complete orbital cycle. He 
described the variability in mi in terms of the presence 
of "bumps" that deform the shape of the photospheric 
absorption line profile from that which it would have if 
its surface were rotating as a rigid body. The individual 
"bumps" were found to travel across the line-profile, from 
blue to red, at rates of ~12-17 km/s per hour. The 
analysis of the pattern of moving bumps seen in the Silll 
lines at 4552, 4567, 4574 A led him to conclude that 
they could be accounted for by four nonradial pulsation 
modes. Two of the modes had periods of ^1/12P and 
1/2P, the latter being associated to "the spectroscopic 
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Fig. 4. — Top: Radial velocity measurements of Si III 4552.62 
(open squares) and He II 4685.7 (shifted by +22 km/s; filled-in 
triangles) for the primary component of Spica plotted as a function 
of orbital phase measured from periastron. The curves correspond 
to the predicted RV velocities according to the Aufdcnbcrg ct al. 
(2009) parameters and adopting u) =255° (solid), w =257° (dash 
)and u) =248° (dots). Middle: predicted difference between the 
centroid of a tidally-perturbed and the unperturbed line profiles. 
Bottom: kurtosis of the theoretical line-profiles. The kurtosis of 
the observed line profiles (triangles) measured at an intensity level 
that avoids the contribution from m,2 displays a similar trend. 



equivalent of the ellipsoidal light variability" which is 
due to the tidal distortion. 

Our limited data set is unable to fully address the 
short-timescale variability since we only have a single 
spectroscopic observation set of 8 exposures taken over 
less than 5 minutes on most nights. The two nights with 
16 exposures (27th and 28th) only cover 10 minutes. 
However, we do have two pairs of nights (March 27 th 
& 23 rd , 25" 1 & 21 s *) where there are small orbital phase 
differences between observations corresponding to a few 
hours of stellar rotation. These pairs allow the charac- 
terization of changes in the "bump" structure over short 
time-scales, though with one complete orbit actually in 
between the observations. The night-to-night variability 
can be easily addressed with such high quality data. 

Figure [5] illustrates the Si III 4552. 62 A line in the 9 dif- 
ferent orbital phases covered by our observations, plotted 
on a velocity scale after applying the correction for or- 
bital motion as listed in Table [31 The variability in this 
line is typical of that seen in all other weak lines. Note 
that the absorption from ni2 can be clearly followed in 
these spectra. In the first five profiles (</> =0.1-0.4, from 
bottom to top), the absorption line originating in the 
secondary component, ni2, is observed to move towards 
maximum negative velocity and then return towards the 
rest-frame velocity. In the last two profiles ((f> =0.8472, 
0.8838), m,2 , s contribution may be seen near its maxi- 
mum positive velocities. In the two profiles at <j) =0.652 



0.92 



Si III 4552.62 
1 r 




-300 -200 



-100 100 
RV (km s- 1 ) 



200 300 



Fig. 5. — Montage of the observed Si III line profiles on a velocity 
scale centered on 4552.62 A, and corrected for the orbital motion 
of the primary using RV values listed Table [3] The spectra are ver- 
tically arranged by orbital phase. The tick marks show the troughs 
of traveling bumps and the secondary's component is marked with 
m 2 . 



the contributions from mi and m-i overlap. 

We have marked in Figure[5]the location of the discrete 
absorptions that are resolved in our data and whose ra- 
dial velocity variations could be measured. In general, 
at least 5 such troughs may be identified, although at 
times some of these may blend together leading to the 
appearance of only 3 troughs (for example, at orbital 
phases 0.3623, 0.3776 and 0.3987). If we assume that the 
troughs retain their identity between two consecutive or- 
bital phases 2 , then we find that between <\> =0.1167 and 
0.1387, the rate by which they travel from blue to red is 
AVbump/ At =22 km/s per hour; and between cf> =0.8472 
and 0.8838, AV bump /At =8 km/s per hour. For <j> =0.35- 
0.39 one of the bumps appears to remain stationary, the 
second moves slightly "red- ward", and the third bump 
seems to split into at least two absorption- like features. 

Thus, the traveling bumps described by Smith (1985b) 
are clearly visible in Figure [5] However, only the troughs 
for <fi ~0.1 in Figure [5] travel at an average displacement 
speed similar to the ~17 km/s per hour speed that may 
be inferred from Smith (1984a; Figures 1-4). At the 
other orbital phases covered by our data the velocity is 
significantly slower. These speeds, of course, are the re- 
sult of the actual motions on the stellar surface combined 
with the projection of the stellar rotation velocity along 
the line-of-sight to the observer. Further insight into the 
phenomenon requires the use of a theoretical model. 

The stability on ~2 hour timescales may be fully ap- 
preciated in Figure [5] In three of the four comparisons 
of Figure [5] there is excellent agreement and the identity 

2 Although it is important to keep in mind that the two con- 
secutive orbital phases in our data correspond to different orbital 
cycles 
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Fig. 6. — Comparison of observed Si III 4552A line profiles, illustrating the changes that occur over A(f> ~0. 02-0. 03 (2-3 hours). Each 
panel shows two observations at similar phases, though separated by several days. The velocity scale is shifted for each profile according to 
the RV values listed in Table [3] The general shape of the profile is very consistently reproduced after one or two full orbits but the small 
scale structure of the bumps is sometimes remarkably consistent (top right) and sometimes highly variable (bottom left). The dates and 
orbital phases for each observation are noted in each panel with the solid curve being the first date/phase. 



of individual bumps is almost entirely preserved. This is 
despite the fact that the observations are separated by 
a full orbital cycle and that diverse variability sources 
could possibly be present. The most striking difference 
appears in the third comparison, corresponding to the 
line profiles from March 21st and 25th obtained at nearly 
the same orbital phase, 0.6515 and 0.6523. These obser- 
vations have the closest phase overlap of all four shown 
in the Figure and yet they present the most significant 
line profile differences. This indicates that the variability 
is not fully orbital-phase locked. 

Figure [5] also provides a view of the manner in which 
the absorption line wings rise from the core to the con- 



tinuum level. There is a marked difference between the 
shape of the line-wings at different orbital phases. For 
example, at orbital phases ~0.3-0.4, the red wing gives 
the profile a more "boxy" shape than at other phases, 
changes that may be quantified using the kurtosis, which 
describes the degree to which a distribution is centrally 
peaked or flat. The bottom panel of Figure [4] shows 
the phase-dependent behavior of the kurtosis, confirm- 
ing the difference that is apparent from Figure [HI in the 
line shapes. The kurtosis was computed over a wave- 
length range that avoids the contribution from mi, and 
the error bars correspond to an uncertainty ±0.2 A on 
the definition of each line wing. 
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TABLE 4 

Summary of model parameters 



Model 


ml 


m2 


Rl 


i 


Po 


V 


e 




dR/Rl 


All sets 10.25-11.8 


6.1-6.97 6.84-8.1 54- 


-66 1.29- 


-2.07 0.009 


-0.058 0.016- 


-0.146 209- 


-258 0.01-0.087 


C 


10.90 


6.80 


8.10 


66 


1.30 


0.028 


0.146 


248 


0.05 


D 


10.90 


6.80 


8.10 


66 


1.30 


0.028 


0.146 


248 


0.06 


E 


10.90 


6.80 


8.10 


66 


1.30 


0.028 


0.146 


248 


0.07 


G 


10.90 


6.80 


8.10 


66 


1.30 


0.018 


0.146 


248 


0.05 


GGD 


10.25 


6.97 


6.84 


54 


2.07 


0.028 


0.067 


255 


0.07 


GGD3 


10.25 


6.97 


6.84 


54 


2.07 


0.028 


0.067 


209 


0.07 


GGD4 


10.25 


6.97 


6.84 


54 


2.07 


0.018 


0.067 


255 


0.07 


GGD6 


10.25 


6.97 


6.84 


54 


2.07 


0.018 


0.067 


255 


0.07 


GG2 


10.90 


6.80 


8.10 


66 


1.30 


0.018 


0.146 


258 


0.07 


H 


10.90 


6.80 


8.10 


66 


1.30 


0.010 


0.146 


248 


0.07 



For the models in this Table, the polytropic index was 1.5. However, we have run models with an index of 1.5, 1.67 and 3. 



In summary: Two categories of line-profile variability 
are present in our data: 1) there are traveling bumps 
that retain their identity over timescales of at least a 
day, though some may change significantly; and 2) there 
is a varying shape of the line-wings which go from being 
very extended to "box-like" on a nightly timescale. We 
will show below that both of these characteristics are a 
consequence of surface flows that can be produced by the 
tidal interaction between the two components. 

4. THE TIDAL FLOW MODEL 

A system is in synchronous rotation when the orbital 
angular velocity Q equals the angular velocity of rota- 
tion cj . In eccentric orbits, the degree of synchronicity 
varies with orbital phase. We use periastron passage as 
the reference point for defining the synchronicity param- 
eter, flper 3 - When at any orbital phase (3 ^1, non-radial 
oscillations are excited, driven by the tidal interactions. 
We refer to the azimuthal components of the forced os- 
cillations as "tidal flows" . 

In the following paragraphs we summarize the general 
characteristics of a simple model developed to gain some 
insight into the effects that the tidal flows introduce in 
the photospheric line profiles. A more detailed descrip- 
tion is provided in Moreno & Koenigsberger (1999) and 
Moreno et al. (2005). 

Consider a binary system in which the light of the pri- 
mary star with mass m\ and radius R\ dominates the 
shape and intensity of a photospheric spectral line. The 
main body of the star, below the surface layer of thick- 
ness Ai? is assumed to have a constant spin angular ve- 
locity u>o with respect to an inertial frame. The binary 
companion has mass and an instantaneous orbital an- 
gular velocity £1 in its (in general) elliptic orbit around 
m\. The axis of stellar rotation is assumed to be per- 
pendicular to the orbital plane. We define (3q — uio/tto, 
where fio is the orbital angular velocity at periastron. 
flo specifies whether the stellar rotation and the orbital 

3 /3 P er = uio/O per = 0.02 P ^jf^^~yTj2 The orbital eccentricity 

is denoted e, the rotation velocity v ro t is given in km/s, the orbital 
period P is given in days, and the stellar equilibrium radius Ri is 
given in solar units. 



periods are synchronized (/?o = 1) or whether the rota- 
tion is super- or sub-synchronous (/3o > 1 or (3q < 1, 
respectively). If /3o ^1.000, then tidal flows on the stel- 
lar surface are induced and non-radial oscillations are 
excited. Our model is based on the fact that the surface 
layer is the one most strongly affected by the external 
gravitational potential of the companion star (see, for 
example, Dolginov & Smel'Chakova 1992; Eggleton ct 
al. 1998) and we therefore compute the response of this 
layer to all the forces in the system, while the rest of 
the star below it is treated as a rigidly-rotating body. 
We find a posteriori, that the dominant contribution to 
the photospheric line-profile variability is produced by 
the azimuthal component of the velocity field and not 
the radial component. In our model, the azimuthal flow 
structure is governed by the tidal forcing and appears 
to be very weakly dependent on the radial component of 
the velocity field. 

A non-uniform Lagrangian grid is constructed over the 
stellar surface by specifying the number of surface ele- 
ments along the equatorial belt and the number of lat- 
itudes for the computation. The number of surface ele- 
ments for each latitude is defined so that in the unper- 
turbed star all elements have similar sizes in the longi- 
tudinal coordinate. For each surface element in the grid, 
an equation of motion is specified in which the follow- 
ing forces are included: the gravitational field of mi and 
rri2, gas pressure, Coriolis force, centrifugal force and the 
viscous shear forces associated with the interface of the 
thin surface layer with the inner stellar body and the 
interface between two adjacent surface elements. The 
surface layer may be treated either in the polytropic or 
in the isothermal approximation. The cases described 
below correspond to calculations with polytropic indices 
of 1.5. Cases run with a different polytropic index or in 
the isothermal approximation yield similar global char- 
acteristics. The simultaneous solution of the equations 
of motion for all surface elements yields values of the 
radial, v r , and azimuthal, v v , velocity fields over the stel- 
lar surface. Azimuth angle, if, is measured from the line 
joining the centers of mi and mi, in the non- inertial 
reference frame centered on mi, and rotating with the 
instantaneous orbital angular velocity £1. Polar angle, 9, 
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is measured from the pole of mi to its equator. 

The model is not restricted to slow rotation nor to 
small excentricities. The only a priori assumptions built 
into the calculation are: 1) only the thin outer layer re- 
sponds to the tidal perturbation; 2) the stellar rotation 
axis is perpendicular to the orbital plane; and 3) the re- 
sponse of the outer layer is either polytropic or isother- 
mal. It must be kept in mind, however that because of 
the single-layer approximation, we compute only the mo- 
tion of the center of mass of each surface element, thus 
ignoring the details of the inner motions within these 
elements, and thus the buoyancy force is neglected. 

4.1. Input parameters 

The input parameters that determine the physical 
characteristics of the tidal interaction effects are: the 
stellar masses, mi, ni2, the orbital period, P or t>, the ec- 
centricity of the orbit, e, the primary star radius, Ri, 
the synchronicity parameter, /?o, the relative depth of 
the surface layer, AR/Ri, and the kinematical viscosity 
of the stellar material, v. In addition, in order to com- 
pute the projection of the surface velocity field along the 
line of sight to the observer, the longitude of periastron, 
w, and the inclination of the orbital plane, i, are required. 
Finally, the properties of the spectral lines that need to 
be specified are the strength of the absorption line ao 
and the line-broadening parameter, k. The latter is as- 
sociated with the microturbulence velocity and for the 
calculations presented here, we use values in the range 
10-15 km/s. 

We performed model calculations for both sets of pa- 
rameters listed in Tabled as well as for intermediate val- 
ues. The range in parameters that was explored is listed 
in the first line of Table |U However, since the Aufden- 
berg et al. (2009) parameter set yields a closer match to 
our observational RVs (Figure [4j, they are adopted for 
the remainder of this paper. The two free parameters 
that remain to be discussed are the depth of the surface 
layer and the kinematical viscosity. The constraint im- 
posed by the basic assumption of a thin outer layer is here 
adopted to mean that the depth of the oscillating layer, 
AR/Ri <0.1. We computed models with AR/Ri =0.05- 
0.08 to determine the dependence of the results on this 
parameter. 

The kinematic viscosity, v determines the efficiency of 
angular momentum transfer between the surface layer 
and the rigidly-rotating inner body of the star, and it 
also couples the radial motion of each surface element 
with that of its neighboring elements. In order to model 
many of the interesting binary systems, we generally re- 
quire a value of v in the range 10 13 -10 16 cm 2 s _1 to limit 
the velocity amplitudes. When the amplitudes grow ex- 
cessively, neighboring surface elements may become de- 
tached, or they partially overlap, causing the computa- 
tion to halt. Our choice of v is thus determined by the 
smallest value that allows the computation to proceed 
for as many orbital cycles as initially specified. The need 
for a large v may be related to the fact that the oscil- 
lations are restricted to a single surface layer and that 
there is no built-in mechanism for energy to be trans- 
ported out of the layer. On the other hand, high values 
of the viscosity are not unusual in astrophysical plasma 
calculations. For example, for a so-called alpha-viscosity 
disk such that v = acH, where c is the sound speed and 



H is the disk scale height, using a typical value a =0.01 
for T Tauri disks (Hartmann et al. 1998) and a typical 
disk structure (D'Alessio et al. 1999) one derives viscosi- 
ties precisely in the same range 10 13 -10 16 cm 2 s _1 . 

0.15 r, . . . . 1 . . . . 1 . . . , 1 — . — , — r. 
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Fig. 7. — The azimuthal velocity field obtained from Scharle- 
mann's (1981) analytical formulation (continuous curves) for a bi- 
nary system with mi = 10.9 Mq, ni2 = 6.8 Mq, Ri = 8.1 Rq, 
P = 1.001 for latitudes from bottom to top 0° (equator), 9°, 17°, 
26°, 35°, 44°, 52°, 61°, 70°, 93°. The results for the same binary 
parameters obtained from our one-layer calculation are illustrated 
with the dash curves. 



One final remark concerns the synchronization param- 
eter /3o- In an eccentric orbit, (3q changes its value as 
a function of orbital phase, being smallest at periastron 
(where the orbital motion is fastest) and largest at apas- 
tron. The values of (3q that we quote to identify a par- 
ticular model calculation refer to its value at periastron 
passage, although its changes over the orbital cycle are 
included in the numerical calculations. From an obser- 
vational standpoint, it is not straightforward to deter- 
mine Po since it requires knowledge of the rotation rate 
of the star in the layers underlying the visible portion of 
the stellar atmosphere. As will be shown below, the ef- 
fects produced by the tidal interaction modify the shape 
of the photosphcric absorption-line profiles. Hence, al- 
though the rotational velocity of the star is in principle a 
quantity that is derived from observations, it is not clear 
whether the value derived from the line profiles corre- 
sponds to the actual rotation speed of the inner stellar 
regions. The uncertainty in vsini translates into an un- 
certainty in j3o- 

As the code is a time-marching algorithm, we let the 
simulations run for at least ~45000 time steps (roughly 
corresponding to 30 orbital cycles), long enough to get 
rid of the perturbations generated by the initial condi- 
tions. Once these initial perturbations disappear, the 
solution is stable over many more orbital cycles for all 
model parameters shown here. 

4.2. The azimuthal velocity field 

We use the term "tidal flows" to refer to the large-scale 
azimuthal velocity patterns that are induced by the grav- 
itational field of the binary companion. An intuitively 
simple way to grasp the essence of the idea of tidal flows is 
given by the description of Tassoul (1987): "Assume that 
at some initial instant the primary is rotating rigidly with 
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a constant supersynchronous angular velocity. If there 
were no companion, the system would remain forever 
axisymmetric. Because of the presence of the secondary, 
however, a fluid particle moving on the free surface of the 
primary will experience small accelerations and decelera- 
tions in the azimuthal direction." In the simplest of cases, 
the tidal flow pattern may be characterized in terms of 
four large contiguous zones over which the acceleration 
alternates between positive and negative values. This de- 
fines four quadrants over the stellar surface, the first of 
which always lies near the sub-binary longitude, and in 
which the acceleration is negative for supersynchronous 
rotation. Figure [7] is an example of such velocity pertur- 
bations obtained from the analytical solution given by 
Scharlemann (1981) for a binary system with mi =10.9 
Mq, m 2 =6.8 M Q , R x =8.1 R Q , and P=4.01d. The Fig- 
ure shows a plot of the azimuthal velocity vs. azimuth 
angle for different latitudes, starting near the equator. 
Because his formulation is valid only for circular orbits, 
near-synchronous rotation and inviscid fluids, we used 
=1.001 and v =0.04 Rg/day, the minimum allowed 
viscosity by the code for this case. The dashed curves in 
this Figure illustrate the results from the calculations of 
our code for a AR/Ri =0.025 surface layer and the same 
orbital and stellar parameters. The one-layer velocity 
perturbations are very similar in shape and in amplitude 
to those of Scharleman's model at latitudes lower than 
~44°. At intermediate latitudes (~44°-70°) the shape 
is still the same, although the one-layer amplitudes are 
systematically larger. 
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Fig. 8. — An example of the azimuthal velocity perturbations at 
the equator from one of our eccentric binary models (CaseGGD6). 
These azimuthal perturbations do show the typical two-lobe tidal 
morphology but are significantly more complex than those shown 
in the previous Figure due to the larger eccentricity and rotation 
rates. The usual two-lobe tidal structure is always reproduced but 
with significant smaller-scale structure induced by the non-zero 
eccentricity and larger stellar rotation rates. 



The velocity field becomes significantly more complex 
when eccentricity and faster stellar rotation rates are in- 
troduced, as illustrated in Figure [H There we show the 
azimuthal velocity field for the GGD6 one-layer model 
calculations for one of the parameter sets similar to those 
of Spica. Note that the pattern still consists of two max- 
ima and minima, but their shape can no longer be de- 
scribed as simply sinusoidal and, in addition, a great 



deal of substructure is apparent. Furthermore, because 
the orbit is eccentric, the pattern changes from one or- 
bital phase to another. The maxima appear around az- 
imuth angle ^90-150° and ~270-330°, and correspond 
to material that is flowing in the direction of the stel- 
lar rotation. The minima correspond to surface material 
that, with respect to the rigidly-rotating inner region, is 
flowing in the opposite direction; with respect to the ex- 
ternal observer, this material is flowing slower than its 
corresponding rotation speed. 

These "tidal flows" are represented for our model run 
GGD in Figure [§] (top left) with light and dark col- 
ors, corresponding to positive and negative perturbation 
speeds, respectively. This Figure also includes maps of 
the radial perturbations in both velocity, v r (middle), 
and in the stellar radius, R v (right). 

It is important to note that the "tidal bulge" is usu- 
ally characterized as a region over which the maximum 
radial extent of the star is attained. Our model suggests 
that, although on average the tidal bulge region is more 
extended, the detailed structure of the tidal bulge for 
an eccentric, super-synchronously rotating star such as 
Spica consists of small regions of alternating large and 
small radial extent. These are accompanied by similar 
oscillations in the velocity perturbations. Hence, we con- 
clude that even for a relatively small eccentricity, the 
"dynamical tide" plays a major role and cannot be ne- 
glected. In Figure fTUl we illustrate the azimuthal velocity 
field evolving over an orbit on the star's hemisphere that 
would be visible to a distant observer. The velocity fields 
are shown for four different orbital phases: 4> =0 (peri- 
astron), 0.25, 0.5 (apastron) and 0.75, and assuming a 
longitude of periastron of uj =255°. 

A second important point to note is that the pattern 
of flow velocities is relatively fixed with respect to the 
coordinate axis that rotates with the binary system; it is 
not fixed with respect to the coordinate axis that is ro- 
tating with mi. Thus, this is conceptually very different 
from the idea of "spots" that have relatively fixed loca- 
tions on the stellar surface and rotate in and out of our 
field of view with the rotation angular velocity. Rather, 
the actual pattern changes with the orbital period; but 
given that there are two maxima and two minima per 
cycle, this leads to variability on timescales 1/2 P. 

4.3. Line profile variability 

Once the surface velocity perturbations have been 
computed, the projection of the velocity field along 
the line of sight to the observer provides the Doppler 
shift needed to produce the rotationally-broadened and 
tidally-perturbed line-profile that may be compared with 
observations. The equations needed to perform the pro- 
jection along the line-of-sight to the observer are de- 
scribed in Moreno et al. (2005). The Doppler shifts 
derived from these projections are then applied to the 
"local" line profiles; i.e., the shape of the spectral line 
that is assumed to correspond to the emergent line radi- 
ation from each surface element. 

Ideally, the local line profile should be derived from 
a model atmospheres computation. In our current cal- 
culation, we adopt a Gaussian local line profile for the 
following reasons. 

First, we have investigated the effect of various local 
line profile shapes on the shape of the total perturbed 
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Fig. 9. — Azimuthal (left) and radial (middle) velocity perturbations, and radial extent of the surface (right). Maps are color-coded 
such that white corresponds to maximum perturbation and black to minimum. For Av v , maximum (positive) is in the direction of stellar 
rotation; for v r and R v , maximum correponds to expansion. The accompanying plots provide the values for the equatorial belt and are 
centered the same as the maps are at <p =270° . The displayed model is GGD at orbital phase 0.11. 



integrated profile and found the variation to be a very 
small. We used CMFGEN model atmosphere profiles 
(Hillicr & Miller 1998) for T e// =25000K and log(g)=3.75 
(kindly provided by John Hillicr) as well as Gaussian 
and Voigt profile fits to these CMFGEN profiles. The 
CMFGEN model provides emergent flux (local line pro- 
files) for a set of locations across the stellar disk (/z). 
These model atmosphere profiles are then fit with Gaus- 
sian and Voigt functions at every /x so an accurate in- 
tegrated line profile can be made using different local 
line profile shapes at every fx. For the Voigt function 
fits, both parameters arc varied to minimize the devia- 
tion between Voigt and CMFGEN profiles. Typical val- 
ues for the damping parameter were 0.27 to 0.48. This 
illustrates how the variation in local line profile shape 
influences the total perturbed line profile and is detailed 
in the appendix. All local profiles give essentially identi- 
cal rotationally-broadened line widths. The width of the 
total line profile is rotational- velocity dominated and the 
intrinsic line width and shape plays little role. The lo- 
cal line profile shape is somewhat more significant when 
computing the detailed shape of the discrete absorption- 
like features. A wider local profile will tend to smooth 
out these bumps, increasing their width and decreasing 
their amplitude, but the bump velocity and number of 
identifiable bumps remains unchanged for all the cases 
we examined. It is interesting to note that the CMF- 
GEN model spectra produce slightly narrower local line 



profiles for the Si III line (smaller wings, shallower core) 
than do the Gaussian or Voigt profiles and would slightly 
increase the relative contrast of the bumps. We have in- 
vestigated this effect using our computed tidal velocity 
fields as well as simulated velocity fields. These local 
line profile effects were found to be very small and can 
be easily neglected for the velocity fields and perturbed 
line profiles calculated here. See the Appendix for more 
details. 

Secondly, it is important to note that given the dynam- 
ical conditions on the stellar surface, non-uniform tem- 
perature and microthermal speed distributions are most 
likely present. In addition, the shearing flows may lead 
to radial temperature gradients that differ from point to 
point over the stellar surface. Without a full 3D radia- 
tive transfer tratment, the impact of these effects on the 
shape of the local line profiles is difficult to predict and 
constitutes a source of systematic error. Thus, adopting 
local line profiles from a 1-dimensional model atmosphere 
is subject to a wide range of systematic errors and is as 
much of an approximation as is the use of a Gaussian 
local line profile. 

Table [4] summarizes the range in input parameters that 
was explored in the nearly 200 model calculations run 
for Spica as well as some specific examples. A general 
feature of all the calculations is that the line-profile un- 
dergoes a strong change in its shape on a day-to-day 
timescale, transitioning from a Gaussian-like profile to a 
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Fig. 10. — The azimuthal velocity field from the GGD model calculation at four different orbital phases, oriented so as to illustrate the 
hemisphere that would be visible to an observer at these particular times; clockwise, starting at the top left map: phase (periastron), 
0.25, 0.50 (apastron), and 0.75. The rotation axis is tilted by 60° with respect to the line-of-sight to the observer. The white broken 
line corresponds to the sub-binary longitude. White(black) corresponds to maximum azimuthal perturbation in the direction of rotation 
(opposite to rotation). 



more "boxy" shape. These changes may be associated 
with the "equilibrium tide" component of the gravitati- 
nal perturbation. As pointed out by Smith (1985a), it is 
the spectroscopic analogue of the ellipsoidal variability in 
the light curves and produces a period of 1/2 the orbital 
period. 

The change in the line-prohles associated with this ef- 
fect can be quantified using the kurtosis. The bottom 
panel in Figure 3] shows the variation of the kurtosis in 
the computed line profiles for the model Case GGD over 



the orbital cycle. As expected from the above discussion, 
there are two maxima and two minima in this curve, cor- 
responding to the transitions into and out of the "boxy" 
line-profile shape. In order to compare the theoretical 
kurtosis with the observational data, we also measured 
the kurtosis in the Si III 4552. 62A line of our 9 spec- 
tra. Here, however, the wavelength interval measured 
was chosen to avoid the contribution from the binary 
companion. The trends in these measurements are plot- 
ted in Figure 2] as the filled- in trianges. The error bars 
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Fig. 11. — Comparison of observations (dots) with theoretical 
profiles for orbital phases when the m.2 line profile is clearly sep- 
arated from that of mi. The size of the dots corresponds to the 
uncertainty in the data. The theoretical profiles reproduce the gen- 
eral shape of the line- wings as well as the approximate number and 
location of the "bumps" . 



in the data correspond to the uncertainty introduced by 
increasing or decreasing the wavelength interval by O.lA. 
Thus, we conclude that our model adequately reproduces 
the variability trends in the observational line profiles 
showing that tidal flows are an important source of non- 
radial velocities. 

It is interesting to note that the line-profile variability 
introduces a systematic error in the measurement of the 
line centroid, as illustrated in the middle panel of Fig- 
ured! There we plot the difference between the centroid 
measured on the perturbed line-profile and the unper- 
turbed profile at the same orbital phase. For the cases 
considered here, the largest error occurs in the orbital 
phase interval 0.6-0.7. 

A comparison of the observations with examples of the 
predicted line-profiles is shown in Figure ITT1 illustrating 
that the model reproduces the orbital-phase dependent 
variations in both the general shape and in the number 
and approximate location of the discrete absorptions. It 
also reproduces the "blue-to-red" motion of the discrete 
features as illustrated in the montage of all the theoreti- 
cal line-profiles illustrated in the traveling bumps section 
of the appendix. Also interesting to note is the pres- 
ence in our computed profiles of the "red spike" found by 
Smith (1985b) and attributed to the horizontal motions. 
Finally, we note that since the profiles for all orbital 
phases are computed sequentially during the same run, 
the comparison of their time-dependent changes (from 
one orbital phase to the next) with that of the data con- 
stitutes yet another constraint. FigurefTTI illustrates that 
the model calculation is consistent with the data also in 
this respect. 

The one important inconsistency between the model 
line-profile calculation shown in Figure [IT] and the data 
is that we used lo =209° for the imput parameter corre- 
sponding to the argument of periastron, instead of the 



value u> =255° implied by the RV curve. We wish to 
stress, however, that the theoretical profiles presented 
in Figure [11] arc not the product of a fitting process; 
these example profiles result directly from the tidal in- 
teraction computation using one of the sets of stellar and 
orbital parameters that have been derived for the Spica 
system. Small modifications in some of the input pa- 
rameters can have a significant impact on the details of 
the line-profiles. A more in-depth discussion of the ef- 
fects on the line-profiles due to different input parame- 
ters will be presented elsewhere. However, it is impor- 
tant here to illustrate the effects produced by the two 
free parameters of the model, v and dR/R\. Figure fT2l 
shows that the large-scale orbital phase-dependent vari- 
ability is generally much stronger than the variations 
produced by changing v and dR/Ri. Hence, once the 
stellar and the binary parameters are constrained, the 
number of free parameters that enter into our calcula- 
tion are greatly reduced. Finally, Figure [13] illustrates 
the perturbed line-profiles computed with the full veloc- 
ity field (as in Figures [11] and [T2"|) compared with the 
profiles computed with only the horizontal velocity com- 
ponent, Av v (dashes). These profiles are nearly identi- 
cal, showing that the radial component of the velocity 
field contributes very little towards the total line-profile 
variability. 

5. CONCLUSIONS AND FUTURE DIRECTIONS 

In this paper we have presented very high quality ob- 
servations of line-profile variability. The combination of 
high spectral resolution with very high precision allows a 
much more detailed comparison of observation and the- 
ory Though we do not have a continuous sequence of 
spectra, we presented some observations that show re- 
markable similarities after a full orbit while others show 
significant differences between orbits. 

We find that the Moreno et al. (2005) approach to 
computing line profile variability from tidal flows yields 
theoretical profiles that can be very similar to observa- 
tional data. Specifically, the theoretical profiles display 
two characteristics of Spica's line-profile variability: 1) 
the "boxy" appearance of the line that is observed twice 
per orbital cycle and that is associated with the "equilib- 
rium tide" ; and 2) the superposed narrow absorption- like 
features that travel from the blue wing to the red wing 
on timescales of ~1 day and that may be associated with 
the "dynamical tides" . This is also the 'spike' reported 
in Smith 1985b. Both of these characteristics are caused 
primarily by the horizontal components of the surface ve- 
locity field and thus, we conclude that tidal flows are the 
dominant contributor to Spica's line-profile. Because the 
simple one-layer model cannot be expected to yield infor- 
mation on the internal structure of the star, the relatively 
good coincidence between the ab initio calculation and 
the observations indicates that in a close binary system 
such as Spica, the tidal forcing by the companion is the 
dominant factor determining the surface perturbations. 
The azimuthal velocity fields that are calculated from 
tidal interactions have the potential to explain much of 
the observed variability. 

Given the fact that our computation is performed from 
first principles with no built-in assumptions regarding 
the surface velocity field and heavy physical constraints 
on our input parameters, we consider the coincidence 
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Fig. 12. — Computed line profiles showing the effects of changing the depth of the oscillating layer, dR/R\ and the viscosity, v. Model 
runs C, D, E, G, GG2 and H were used to illustrate variation with layer depth and viscosity at two orbital phases - the left panel is <f> =0.09, 
the right is rp =0.85. The dotted lines correspond to the unperturbed, rotationally-broadened line profile and the solid lines are perturbed. 
The layer depths and viscosities for each line are shown in each panel. The top group of line profiles have a fixed layer-depth with a varying 
viscosity while the bottom group varies layer depth. The general shape of the line is most sensitive to the depth of the oscillating layer, 
but the phase-dependent effects dominate over the relatively smaller differences introduced by the choice of this parameter. 
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Fig. 13. — Line profiles for the Spica GGD model computed with 
the total velocity field compared with the profiles computed with 
only the horizontal velocities, Av v , showing that the radial com- 
ponent contributes very little to the line-profile variability. Shown 
are orbital phases as in Figure [5] 



between the theoretical and observational line profiles to 
be very encouraging. We do use a single layer fluid ap- 
proximation but all the dynamics are controlled from the 
fundamental equations of motion. The structure of the 
fluid flow is not constrained to fit any specific mathe- 
matical form. The velocity field, subject to the single- 
layer approximation, is computed directly from the equa- 



tions of motion. This illustrates the importance of prop- 
erly treating surface flows (the horizontal component of 
the surface motions) when analyzing observed line profile 
variability. 

There are several aspects of this model that can be 
readily explored. The strength of the narrow spectral 
features (bumps) are not ideally reproduced. The feed- 
back from shear energy dissipation into surface tempera- 
ture and local line profile shape should be included to test 
the significance of dissipation in computing the predicted 
profiles. Ideally, model atmosphere line profiles and limb 
darkening should be used in the general scheme, although 
a proper treatment requires the use of a 3D model atmo- 
sphere calculation. In the meantime, however, we have 
shown the shape of the local line profiles to be of little 
significance for this particular system. Also, incorporat- 
ing a more realistic stellar structure into the calculation 
should allow exploration of the effects produced by the 
interplay between the tidal forcing and the normal modes 
of oscillation. 
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Fig. 14. — The pre-processing properties of ESPaDOnS data as reduced by Libre-Esprit. Box a shows the order overlap as a function 
of wavelength. Box b shows the dispersion as a function of wavelength. Box c shows the corresponding spectral-resolution at single-pixel 
sampling. Box d shows the number of spectral-pixels used in the processing scripts to bin in wavelength to a constant lA per spectral-pixel 
coverage. Order-overlap and non-linear coverage leads to a bi-modal distribution of binned-pixels as well as a decrease with wavelength. 
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7. APPENDIX 

The non-linear and irregular wavelength sampling of 
the ESPaDOnS spectropolarimeter led us to develop rou- 
tines to linearize the output spectra. Libre-Esprit pro- 
duces wavelength-calibrated spectra for every possible 
spectral pixel in the focal plane. As with all cross- 
dispersed echelle spectrographs, the dispersion, order- 
overlap and sampling are all functions of wavelength. An 
additional complication is that the ccd pixels correspond 
to 2.6km/s but the spectra are output by Libre-Esprit at 
1.8km/s resolution, 30% smaller than a pixel. These are 
called ccd and spectral pixels respectively. Figure [Lfa 
shows the order overlap as a function of wavelength for a 
typical ESPaDOnS exposure in our data-set. The order 
overlap increases from in overlap from 20A at 4000A to 
over 60A around 7000A before falling strongly to become 



a more than 40A gap in wavelength coverage by 10000A . 
The dispersion, shown in Figure fTib . runs from 0.02A per 
pixel to 0.06A per pixel in proportion to the wavelength 
increase. This leads to an essentially constant spectral 
resolution, shown in Figure [T4"c as roughly R~170,000. 
As the full-width-half-max of most calibration-lamp lines 
is over 2 pixels, the spectral resolution is around 70,000. 

We desired a monotonic wavelength grid that assigns 
one intensity to one wavelength and combines all pixels 
covering the same wavelength to achieve higher signal- 
to-noise ratio. The processing scripts we wrote calculate 
the number of spectral-pixels to bin to achieve 0.12A 
per pixel coverage for all wavelengths. This is shown in 
Figure fMh . Since the order overlap is significant, in some 
cases 30% of an order, there is a bi-modal distribution 
in the number of pixels binned. For the part of an order 
with overlapping coverage, a linear wavelength coverage 
gives roughly twice the number of pixels binned. Beyond 
roughly 8000A , the order overlap disappears and the 
dispersion reaches 0.05A per pixel resulting in a nearly 
regular 2-pixel binning factor. 

After applying the rebinning routines, the dispersion 
remains constant at 0.12A per pixel as seen in Figure fT5b . 
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Fig. 15. — The post-processing properties of ESPaDOnS data. Box a shows the rebinned dispersion per pixel. Box b shows the resulting 
spectral resolution at single-pixel sampling. Box c shows the signal-to-noise improvement factor calculated as the square-root of the number 
of binned spectral-pixels. A polynomial fit is overplotted. Box d shows the signal-to-noise for each order calculated by the Libre-Esprit 
script (lower curve). The binning process results in an improved signal-to- noise (upper curve). 



This results in a spectral resolution that is wavelength de- 
pendent, going from 30,000 to 80,000 per spectral-pixel 
at one pixel per resolution element. This is shown in Fig- 
ure [TSb. As stated above, the physical spectral resolution 
of the raw output is roughly 70,000 with over 2 spectral 
pixels per resolution clement. This new binned output 
undersamples 2:1 at 4000 A increasing to 1:1 sampling at 
7000A and rising with longer wavelengths. The signal- 
to-noise will be increased roughly as the square-root of 
the number of pixels averaged as all our observations are 
in the high signal-to-noise regime. The signal-to-noise 
improvement factors are shown in Figure [15b for both 
order-overlap and no-overlap. Libre-Esprit outputs an 
estimate of the signal-to-noise for each order as part of 
the reduction routine. An example of this is seen in the 
lower curve of Figure [T5H . After the regularization of the 
wavelength array the signal-to-noise is increased. Using 
the modest improvement with no order-overlap, the re- 
sulting signal-to-noise is shown as the upper curve in Fig- 
ure I15H . Though these routines result in a wavelength- 
dependent resolution (with constant dispersion), the reg- 
ular, linear, monotonic wavelength sampling allows for 
easier processing of large numbers of observations along 



with generally higher signal-to-noise. 

We also performed some experiments on the accuracy 
of the signal-to-noise predictions given by Libre-Esprit. 
In the very high signal-to-noise regime, systematic effects 
caused by flat fielding issues, non-linearity or instrument 
instabilities can become significant. The Libre-Esprit re- 
duction logs output two estimates of the signal-to-noise 
based on the flux. There is an estimate for each "spec- 
tral pixel" and "ccd pixel" . The difference between these 
pixels is the resolution of 2.6km/s per physical ccd pixel 
compared to 1.8km/s per spectral pixel as resampled in 
Libre-Esprit. 

The noise properties of the ESPaDOnS data in the 
high signal-to-noise regime can be best seen in our best 
observation set: March 27th. These properties are mea- 
sured directly on the Libre-Esprit output data without 
applying our wavelength-regularizing scripts. Figure [16] 
shows four separate spectra of Spica for order number 
44 centered at 5150A. There are two complete polariza- 
tion measurements which each consist of two intensity 
spectra output by Libre-Esprit. Each individual spec- 
trum there is calculated from four individual exposures 
each with two orders. Thus, each spectrum in Figure 
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Fig. 16. — The noise properties of the ESPaDOnS data in the highest signal-to-noise observation set: March 27th. Box a shows the 
spectrum of Spica in order 44 (5150A). We have two complete polarization sets each with two spectra from Libre-Esprit. Each individual 
spectrum here is calculated from two polarized orders in 4 separate exposures. The Libre-Esprit package gives signal-to-noise ratio's of 
1350 / 1600 for each spectral / ccd pixel. Box a shows the four individual spectra after a simple linear continuum normalization offset for 
clarity. Box b shows the difference in wavelength solutions in pico-meters and pixel fractions after removing an offset of 1-2 pixels between 
individual solutions. The wavelength dependence of the shifts across an order are negligible. Box c shows the residual noise after subtracting 
a 3-pixel boxcar smooth from each intensity spectrum. This highlights the noise and removes spectral structure. Box d shows the cross 
correlation between each of the four noise spectra and the top noise spectra from box c. The top curve represents an autocorrelation. A 
clear signature from the boxcar smooth is seen. The lower three plots also show a correlation, but spectrally offset by 1 or 2 pixels. The 
wavelengths differ by less than 10% of a pixel but a 1-2 pixel shift between exposures has been recorded by Libre-Esprit. This correlation 
between the noise spectra is likely caused by structure on the CCD. 



\Wk represents an average of 8 individual spectra. The 
Libre-Esprit package gives signal-to-noise ratio for each 
plotted spectrum of roughly 1350 / 1600 per spectral 
and ccd pixel. The wavelength solution is derived inde- 
pendently for each 4-exposure polarization measurement. 
The shift between the wavelength solutions is typically 
0-2 pixels. In Figure flBb . the change in wavelength solu- 
tion is shown to be a small fraction of a pixel across this 
spectral order. The illumination pattern on the CCD 
changes only slightly between individual exposure sets. 

The noise in the spectra can be estimated by subtract- 
ing a boxcar-smoothed spectrum from a raw spectrum. 
As the full-width-half-maximum of spectral features is 
roughly 2.2 pixels, any wavelength variation present at 
the 1-spectral-pixel level is surely noise. FigurefTfib shows 
residual noise from each spectrum after subtracting a 3- 
pixel boxcar smooth. This highlights the noise and re- 
moves spectral structure. This noise-spectrum can then 



be cross-correlated to determine the level of any system- 
atic errors present. Figure I16H shows the cross corre- 
lation between each of the four noise spectra and the 
first noise spectra (the top spectrum in box c). The top 
curve represents an autocorrelation and a clear signa- 
ture from the boxcar smooth is seen. The lower three 
cross-correlation curves also show this same signature, 
but spectrally offset by 1 or 2 pixels and at significantly 
lower amplitude. The wavelength solutions differ by less 
than 10% of a pixel across an order and are generally 
stable. However, a 1-2 pixel shift between exposures has 
been recorded by Libre-Esprit between the four wave- 
length solutions . This correlation between the noise 
spectra is likely caused by structure on the CCD. 

8. TRAVELING BUMPS 

In this appendix section we will describe some general 
properties of the model calculations. Figure [T7l shows a 



a Vir Spectral 



TABLE 5 
ESPaDOnS Orders 



N 
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A c 


A l 


N 


\n 
A 


\ 

A c 




61 


369.12350 


372 


375.10460 


11 


542.71300 


552 


560.72960 


60 


372.96250 


377 


381.41990 


40 


556.02290 


566 


574.94500 


59 


379.20690 


384 


387.95310 


39 


570.00300 


581 


589.90050 


58 


385.66620 


390 


394.71350 


38 


584.70170 


596 


605.65430 


57 


392.34930 


397 


401.71090 


37 


600.17300 


612 


622.27430 


56 


399.26270 


404 


408.96230 


36 


616.49410 


629 


639.82890 


55 


406.42750 


412 


416.47950 


35 


633.71090 


647 


658.39830 


54 


413.85750 


419 


424.27950 


34 


651.92970 


666 


678.07840 


53 


421.55920 


427 


432.37580 


33 


671.44260 


686 


698.79780 


52 


429.55120 


435 


440 78550 


32 


692 43090 


708 


720 63430 


51 


437.85630 


444 


449.52940 


31 


714.77320 


731 


743.87970 


50 


446.48130 


453 


458.62520 


30 


738.60490 


755 


768.67450 


19 


455.45800 


462 


468.09680 


29 


764.07990 


781 


795.17950 


18 


464.79900 


472 


477.97030 


28 


791.37460 


809 


823.57770 


17 


474.53390 


482 


488.26610 


27 


820.69110 


839 


854.07930 


16 


484.68840 


492 


499.01680 


26 


852.26270 


871 


886.92720 


45 


495.27700 


503 


510.24760 


25 


886.36000 


906 


922.40310 


41 


506.33940 


515 


521.99910 


24 


923.29860 


944 


960.83560 


43 


517.91600 


527 


534.30210 


23 


963.44940 


985 


1002.6098 


42 


530.02500 


539 


547.19600 


22 


1007.2513 


1029 


1048.1787 



This Table shows order number, calculated beginning wavelength, 
Libre-Esprit log "central" wavelength and calculated ending wave- 
length. 




Velocity Velocity 

Fig. 17. — This illustrates the time evolution of the line pro- 
file over a single orbit. The left panel shows deviations from the 
average profile over an orbital cycle. The right hand panel shows 
the average line profile as well as some individual line profiles from 
various orbital phases. The velocity is corrected for mi's orbital 
motion. 



typical model run. The code is allowed to run for roughly 
30 orbital cycles to damp out any perturbations from 
initial conditions. The perturbed velocity field is then 
calculated at each timestep, as is its projection along the 
line-of-sight to the observer, from which a total line pro- 
file is output. The right hand panel of Figure [T71 shows 
the average line profile over the orbital cycle as the dark 
solid line. Several individual spectra from throughout 
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the orbit are also over-plotted as the dashed and dotted 
lines. Once the average line profile is subtracted from 
each individual measurement, a grey-scale plot of the 
deviations can be constructed to show how individual 
surface perturbations move across a line profile and pre- 
serve their identity over an orbital cycle. The left hand 
panel is dominated by the main tidal bulge - the light and 
dark patches that are most prominent on the line wings 
(largest velocities). However, the smaller-scale surface 
perturbations can be more easily seen near line center as 
the individual diagonal streaks. 




4553.5 4554.0 4554.5 

Wavelength In Angstroms 



Fig. 18. — This illustrates the local line profiles used in comput- 
ing perturbed average spectral line profiles. The solid line repre- 
sents the CMFGEN local line profiles. The dashed line is a Voigt 
profile fit to the CMFGEN line profile. The dot-dash line is a 
Gaussian fit to the CMFGEN line profile. Fits are shown for both 
disk-center (fJ.=l) and for the limb (fi=0). The difference between 
Voigt and Gaussian fits lie primarily in the width and the wings of 
the profiles. 



9. LOCAL LINE PROFILE EFFECTS 

In order to test the accuracy and utility of using Gaus- 
sian local line profiles, we performed an investigation 
as to the effect of using different shape local-line pro- 
files in computing the integrated rotationally broadened, 
tidally-perturbed profiles. We obtained CMFGEN model 
atmosphere profiles with T e //=25000K and log(g)=3.75 
provided by John Hillier. Since the CMFGEN model 
provides emergent flux (local line profiles) for a set of 
projected locations on the stellar disk (/i) this database 
of profiles represents our starting point to simulate local 
line profile effects. We now have a model local line profile 
shape for any point on the stellar surface. 

These model atmosphere profiles are then fit with both 
Gaussian and Voigt functions at every fi. These fits 
are straightforward using IDL's curve-fitting routines: 
Gaussfit and Curvefit. Once we have CMFGEN, Gaus- 
sian and Voigt local line profiles for every \i an accurate 
integrated line profile can be made to simulate an obser- 
vation of a tidally perturbed star. The projected velocity 
is known for every surface element as is the value for fi 
and emergent flux at every surface element. An example 
of these local line profiles is shown in Figure [TH] 

To illustrate the effects of local line profile shapes on 
the perturbations caused by macroscopic surface flows, 
we computed a perturbed line profile for a velocity field 
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using all three of these local line profile shapes. A simu- 
lated perturbing velocity field is created as a combination 
of spherical harmonic functions and random fluctuations 
with a maximum amplitude of ±15km/s. This creates a 
highly perturbed line profile with known properties that 
can clearly demonstrate the effect of local line profile 
shapes. Figure [TH shows a perturbed velocity field in 
the top panel as well as integrated line profiles in the 
bottom panel. The integrated line profiles are calculated 
using the perturbed stellar velocity field. The CMFGEN 
(solid), Gaussian (dot-dash) and Voigt (dash) local Si 
III line profiles are used to generate an integrated pro- 
files seen in the bottom panel. There is essentially no 
difference between integrated line profiles even though 
significant differences in local line profile shapes are ap- 
parent. 



expected, the difference in local line profile shape was 
negligible. Figure [20] shows integrated line profiles com- 
puted for Case GGD3 parameters. The more boxy local 
line profiles are used for the solid curve and the Gaussian 
local line profiles are used for the dotted curve. There is 
a very small difference in the depth of the perturbations 
but the overall morphology is very similar. 

We consider the use of Gaussian local line profiles to 
be entirely sufficient to demonstrate this tidal-flow model 
for integrated line profiles. The integrated line profile 
intensities differ at the 0.05% level for various local line 
profile shapes. This is considered to be entirely below any 
systematic errors or errors caused by uncertainties in the 
physical parameters used in calculating the properties of 
the binary systems. 
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Fig. 19. — This illustrates the perturbed integrated profiles us- 
ing different local profile shapes. The top panel shows a veloc- 
ity perturbation on the stellar surface. The bottom panel shows 
the integrated line profile when calculated using CMFGEN (solid), 
Gaussian (dot-dash) and Voigt (dashed) local line profiles. The 
difference between local line profile morphology is significant - the 
Gaussian and Voigt local line profiles show broader wings and a 
sharper core while the CMFGEN spectra are more 'boxy'. How- 
ever, these local differences have an entirely negligible effect on the 
integrated profile. 



As a third test of local line profile shape effects on in- 
tegrated line profiles, we adapted our tidal-flow code to 
compute integrated line profiles using a more boxy lo- 

4 2 

cal line profile shape - ~ e~ x instead of ~ e~ x . As 
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Fig. 20. — Computed line profiles with Case GGD3 parameters 
showing that the intrinsic shape of the adopted line profile only 
affects the the depth of the sharp features, but does not modify 
the general characteristics of the perturbed line profile. The dotted 
profile was computed with a Gaussian intrinsic line-shape, while 

4 

the continuous one with a more "boxy" (~ e~ x ) intrinsic line- 
shape. 



10. VELOCITY PERTURBATIONS 

We feel it is worth pointing out that the amplitude and 
size of the 'bumps' fairly tightly constrain the size and 
amplitude of the corresponding velocity fluctuations on 
the stellar surface. In order to produce the bumps we ob- 
serve, one simply needs a velocity perturbation of order 
10-30km/s with a few 'patches' of velocity perturbation 
on the stellar surface for each observed bump. Figure |2"T1 
illustrates this idea. 

We created radially projected velocity perturbations 
on the stellar surface with a random number generator 
scaled to a peak amplitude of ±15km/s. The surface 
of the star was perturbed with this random field scaled 
to have an increasing number of 'granules' on the stel- 
lar surface. This illustrates the effect of decreasing the 
size-scale of velocity perturbations on the integrated line 
profiles. Our experiments were run with spatial scales of 
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Fig. 21. — This illustrates how randomly generated velocity perturbations with various spatial distributions can give rise to line profiles 
similar to those observed. The top panel shows various velocity fields - the 160km/s rotational broadening in the top left and randomly 
generated fluctuations with 2, 5, 10, 20, 50, 100, 200 and 400 fluctuations across the equator. The line profiles corresponding to these 
velocity fields (scaled to ±15km/s) are shown in the bottom panel. Lines with several 'bumps' of roughly the amplitude observed are seen 
in the examples with 10 to 20 equatorial bumps. 



0, 2, 5, 10, 20, 50, 100, 200 and 400 individual perturba- 
tions across an equatorial slice. 

The top panel of Figure [5T] shows first the rotational 
broadening of 160km/s in the top left followed by these 
perturbed velocity fields. Once these perturbations are 
incorporated with the rotational broadening, the corre- 
sponding line profiles are calculated. The bottom panel 
of Figure |2"T1 shows these perturbed line profiles. A fairly 
simple trend can be seen. A small number of perturbed 
surface patches (2-5) results in small bulk-shape devia- 
tions of the line profile. When the number of surface 
perturbation patches reaches 10-20, line profiles resem- 
bling those observed are seen. As the number of sur- 
face perturbations increase, more bumps appear at in- 



Aufdenberg, J. et al., 2009, in prep. 
Baade, D., 1984, A&A, 134, 105. 

Dolginov, A.Z. & Smel'Chakova, E.V. 1992, A&A 257, 783. 
Donati, J.F. et al., 1997, MNRAS, 291, 658. 
D'Alessio, P. et al, 1999, ApJ, 527, 893. 

Eggleton, P.P., Kiseleva, L. G. & Hut, P., 1998, ApJ, 499, 853. 

Fullerton, A.W., Gies, D.R. & Bolton, C.T., 1996, ApJS, 103, 475. 

Gies, D.R. & Kullavanija, A., 1988, ApJ, 326, 813. 

Hartmann, L., et al., 1998, ApJ, 495, 385. 

Herbison-Evans, D. et al., 1971, MNRAS, 151, 161. 

Hillier, D.J. & Miller, D.L., 1998, ApJ, 496, 407. 

Lomb, N.R., 1978, MNRAS, 185, 325. 

Moreno, E. & Koenigsberger, G. 1999, RMA&A, 35, 157. 
Moreno, E., Koenigsberger, G. & Toledano, O. 2005, A&A, 437, 
641. 

Reid A.H.N, et al., 1993, ApJ, 417, 320. 

Riddle, R. L., 2000, Ph.D. thesis, Georgia State Universtiy. 



creasingly smaller amplitudes until (by 400 equatorial 
patches) there's essentially no visually apparent differ- 
ence between perturbed and unperturbed line profiles. 
This example illustrates a simple constraint on the per- 
turbed velocity field coming directly from the observed 
amplitude and number of bumps - relatively large-scale 
surface perturbations (a few to several patches in an 
equatorial slice) with an amplitude that is a fraction of 
the rotational velocity can cause the observed bumps. 
Since the structure of the observed perturbations (the 
line profile shape) is more or less similar from one orbit 
to the next, the corresponding velocity perturbations are 
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